
use "data/ready4analysis.dta",clear
	drop if mi(redz)
	xtset pidp year	

****************************************************************************************
**** FIGURE D1
****************************************************************************************

	recode egp_yrs (min/-1=0), gen(yrs0) // Shift value of non-treated to 0
	xtreg redz i.treat i.treat#c.yrs0 i.treat#c.yrs0#c.yrs0 c.age, fe vce(robust)

		// Year 4
			di (_b[1.treat]+4*_b[1.treat#c.yrs0]+4*4*_b[1.treat#co.yrs0#co.yrs0]) // 0.07
			di (_b[2.treat]+4*_b[2.treat#c.yrs0]+4*4*_b[2.treat#co.yrs0#co.yrs0]) // -0.07
			test (_b[1.treat]+4*_b[1.treat#c.yrs0]+4*4*_b[1.treat#co.yrs0#co.yrs0])= -(_b[2.treat]+4*_b[2.treat#c.yrs0]+4*4*_b[2.treat#co.yrs0#co.yrs0]) // ns. 
		// Year 8
			di (_b[1.treat]+8*_b[1.treat#c.yrs0]+8*8*_b[1.treat#co.yrs0#co.yrs0]) // 0.04
			di (_b[2.treat]+8*_b[2.treat#c.yrs0]+8*8*_b[2.treat#co.yrs0#co.yrs0]) // -0.12
			test (_b[1.treat]+8*_b[1.treat#c.yrs0]+8*8*_b[1.treat#co.yrs0#co.yrs0])= -(_b[2.treat]+8*_b[2.treat#c.yrs0]+8*8*_b[2.treat#co.yrs0#co.yrs0]) // ns. 
					
		* Plotting the marginal effects (Conditional Effect Plot)
		margins, at(treat=(0 1 2) yrs0=(0(1)8)) contrast(atcontrast(r._at) lincom) noatlegend post

		/*
		marginsplot, recast(line) recastci(rline) yline(0, lcolor(black)) ///
		   x(yrs0)                                                ///  
		   plot1opts(lstyle(none)) ci1opts(lstyle(none))                  /// omit graph for never-mobile
		   plot2opts(lpattern(solid) lwidth(thick) lcolor(red))     	  /// estimate for mobile
		   ci2opts(lpattern(dash) lwidth(medthick) lcolor(red))         /// CI for estimate 
		   plot3opts(lpattern(solid) lwidth(thick) lcolor(blue))     	  /// estimate for mobile
		   ci3opts(lpattern(dash) lwidth(medthick) lcolor(blue))         /// CI for estimate 
		   ylabel(-.3(.1).5, grid angle(0) labsize(medium) format(%3.1f)) /// 
		   xlabel(0(1)8, labsize(medium))                                ///
		   xtitle("Years since transition", size(large) margin(0 0 0 2))    ///
		   ytitle("Change in red", size(large))   title(" ")        ///
		   legend(pos(7) ring(0) row(1) order(2 "95%-CI") size(medlarge))    

		  graph export "temp/fig3_egp_combined.png", as(png)  replace
		*/
			
		* Below saves data for R figs			
		qui parmest, saving("temp/figD1.dta",replace) 

		cap frame drop new
		frame create new
		frame change new
				
				use "temp/figD1.dta", clear
					
				matrix varnames = e(at)
				svmat varnames
				
				replace varnames2 = varnames2[_n+1]
				replace varnames3 = varnames3[_n+1]
				replace varnames4 = varnames4[_n+1]
				
				keep if varnames2==1 | varnames3==1 // Only keep treatment effect
				gen treat=1*varnames2+2*varnames3
					label define treat 0 "Non-treated" 1 "Downwardly mobile" 2 "Upwardly mobile", modify
					label values treat treat
				rename varnames4 year	
				keeporder treat year estimate stderr min95 max95

				save "figs/figD1_egp_combined.dta", replace

		frame change default
		frame drop new
	
